//////////////////////////////////////////////////////////////////////
//                                                                  //
//                                                                  //
//              Copyright (C), 2000  SEIKO EPSON Corp.              //
//              ALL RIGHTS RESERVED                                 //
//                                                                  //
//                                                                  //
//   file name : fmod.c                                             //
//                                                                  //
//   Revision history                                               //
//       00/02/14    first release                 M.Igarashi       //
//                                                                  //
//////////////////////////////////////////////////////////////////////

#include <math.h>
#include <smcvals.h>
 
// ALGORITHM
// 1.|dfX|>|dfX2| fmod (dfX,dfX2)= dfX2 * modf( dfX / dfX2,&d)
//   |dfX2|>|dfX| fmod (dfX,dfX2)= dfX
//   |dfX|=|dfX2| fmod (dfX,dfX2)= 0
//

//  the architecture of double floating point
//
//   0 1          11                  31 32                               63 bit
//   -----------------------------------------------------------------------
//  | |   exponent  |                    fraction                           |
//   -----------------------------------------------------------------------
//
//  |               |                   |                                   |
//  |    12 bits          20 bits       |             32 bits               |
//  |            lower word             |            higher word            |
//
//         bit    0         sign bit         (  1 bit  )
//              1 - 11      exponent part    ( 11 bits )
//             12 - 63      fraction part    ( 52 bits )
//
//
 
double fmod (double dfX, double dfX2){

	long lX,lX2;
	unsigned long ulLx2,ulHx2;
	double dfTemp,dfRet,dfMod;
	
	GETW_L(lX,dfX);		// get low
	GETW_L(lX2,dfX2);		// get low
	GETW_H(ulHx2,dfX2);		// get high
	
	ulLx2=lX2&0x7fffffff;	// mask sign

	if(ulLx2==0x00000000 && ulHx2==0x00000000 ){		// dfX2 ==0
	    errno=EDOM;                         // 33: domain error
	    dfRet=NAN._D; 
		return dfRet;
	}

	if(((lX2&0x80000000) !=0 ? -dfX2:dfX2) < ((lX&0x80000000) !=0 ? -dfX:dfX)){		// |dfX|>|dfX2|
		dfTemp=dfX/dfX2;
		dfTemp=modf(dfTemp,&dfMod);
		dfRet=dfX2*dfTemp;
		return dfRet;
	}else{
		if(((lX2&0x80000000) !=0 ? -dfX2:dfX2) > ((lX&0x80000000) !=0 ? -dfX:dfX)){	// |dfX2|>|dfX|
			dfRet=dfX;
			return dfRet;
		}else{				// dfX==dfX2
			dfRet=0.0;
			return dfRet;
		}
	}
}
